function [allsols,  scrs, bestslns] = nonlinearOptimizeselfcalibnormMODREDUCEDTIME(TF,w,h,numParams,initialWeights,numTries,kstart)


numF=size(TF,2);

if (nargin <5)
    initialWeights=ones(numF,1);
end
if (nargin <6)
    numTries=1;
end
if (nargin <7)
   kstart=[w 0 w/2; 0 w h/2; 0 0 1];
end



allsols= repmat({zeros(3)}, 1, numTries);
scrs=zeros(numTries,1);
bestslns=zeros(3,3);
initialSolutions=createInitialSolutions(numTries,w,h,numParams,kstart);
[lb,ub]=findBoundsonK(w,h,numParams);
bestscore=100000000;

for i=1:size(TF,2)
    TF{1,i}=TF{1,i}/norm(TF{1,i});
end

for i=1:numTries
    vi=randperm(size(TF,2));
    if(numParams>2)
       qin=vi(1:2); 
    else
       qin=vi(1:3); 
    end
    [x,error]=nonlincustomFunction(initialSolutions{1,i},{TF{1,qin}},initialWeights,w,h,lb,ub);
    
   allsols{1,i}=convertXTOKselfK(x,w,h);
    % display([ num2str(i) ' out of ' num2str(numTries)]);
    scrs(i,1)=error;
    
    if(scrs(i,1)<bestscore  )
        
        bestscore=scrs(i,1);
        bestslns=allsols{1,i};
        
    end
    
end
end

function [x,error]=nonlincustomFunction(x0,F,weights,w,h,lb,ub)
% according to the discussion weights need to be cholesky decomposed
CholeskyinitialWeights = sqrt(weights);
CholeskyinitialWeights=1000*CholeskyinitialWeights*(1/sum(CholeskyinitialWeights));

levmaroptions=[1E-5, 1E-25, 1E-25, 1E-30, 1E-06]; % this is from the demos
 tolx=1e-14; tolf=1e-14;
%optionsfmincon = optimset('Algorithm','interior-point','Display','iter-detailed','TolFun',tolf,'TolX',tolx);
optionsfmincon = optimset('Algorithm','interior-point','Display','off','TolFun',tolf,'TolX',tolx);
measurements=0;

sizeFs=size(F,2);
numParams=size(x0,2);

type=2;


if(type==1)
    if(numParams>sizeFs)
        measurements= zeros(numParams,1);
        outputformat=2 ;
    else
        outputformat=1  ;
        measurements= zeros(sizeFs,1);
    end

    [ret, x, info]=levmar('computerEssentialErrorSVDNFramesGeneralMOD',  x0, measurements, 2000, levmaroptions,'bc',lb, ub, F,CholeskyinitialWeights,w,h,outputformat);
    info
    ret
    error=sum(abs(info(2)));
elseif(type==2)
    outputformat=5  ; % zzz not sure if i should square each element or the whole thing
    f = @(x)computerEssentialErrorSVDNFramesGeneralMOD(x,F,CholeskyinitialWeights,w,h,outputformat);
  try
 %     x0
    [x,fval,exitflag,output] = fmincon(f,x0,[],[],[],[],lb,ub,[],optionsfmincon);
  %  x
   %  display('_________________________'); 
  catch err
     display('*********** caught some error'); 
     x=x0;
     fval=10;
     save('badcrashfmincon.mat','f','x0','lb','ub','optionsfmincon');
  end
    error=sum(fval);
end


end

function [initialSolutions]=createInitialSolutions(numTries,w,h,numparams,KSTART)

K=KSTART;

fstd=400;
skewmax=2;
initialSolutions=repmat({zeros(1,numparams)}, 1, numTries);


for i=1:numTries
    Knew=K;
    if(i>1)
        Knew(1,1)=KSTART(1,1)+randn()*fstd;% gaussian
        Knew(2,2)=KSTART(2,2)+randn()*fstd;% gaussian
        Knew(1,3)=rand()*w;
        Knew(2,3)=rand()*h;
        Knew(1,2)=rand()*skewmax;
        
        
    end
    initialSolutions{1,i}=convertKTOXselfK(abs(Knew),numparams);
end

end

function [lb,ub]=findBoundsonK(w,h,numParams)
fmax=2000;  fmin=40;
skewmax=10; skewmin=-1;
xmax=w;     xmin=0;
ymax=h;     ymin=0;
armax=1.2;  armin=0.8;
KMAX = [fmax   skewmax        xmax; 0      fmax*armax   ymax;  0            0             1  ];
KMIN = [fmin   skewmin        xmin; 0      fmin*armin   ymin;  0            0             1  ];

lb= convertKTOXselfK(KMIN ,numParams);
ub= convertKTOXselfK( KMAX ,numParams);

end
% /* O: information regarding the minimization. Set to NULL if don't care
%  * info[0]= ||e||_2 at initial p.
%  * info[1-4]=[ ||e||_2, ||J^T e||_inf,  ||Dp||_2, \mu/max[J^T J]_ii ], all computed at estimated p.
%  * info[5]= # iterations,
%  * info[6]=reason for terminating: 1 - stopped by small gradient J^T e
%  *                                 2 - stopped by small Dp
%  *                                 3 - stopped by itmax
%  *                                 4 - singular matrix. Restart from current p with increased \mu
%  *                                 5 - no further error reduction is possible. Restart with increased mu
%  *                                 6 - stopped by small ||e||_2
%  *                                 7 - stopped by invalid (i.e. NaN or Inf) "func" values; a user error
%  * info[7]= # function evaluations
%  * info[8]= # Jacobian evaluations
%  * info[9]= # linear systems solved, i.e. # attempts for reducing error
%  */
